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ABSTRACT 

Aims. Two inter-related peaks occur in high-frequency power spectra of X-ray lightcurves of several black-hole candidates. We further explore 
the idea that a non-linear resonance mechanism, operating in strong-gravity regime, is responsible for these quasi-periodic oscillations (QPOs). 
Methods. By extending the multiple-scales analysis of Rebusco, we construct two-dimensional phase-space sections, which enable us 
to identify different topologies governing the system and to follow evolutionary tracks of the twin peaks. This suggests that the original 
(Abramowicz & Kluzniak) parametric-resonance scheme can be viewed as an ingenuous account of the QPOs model with an internal 
resonance. 

Results. We show an example of internal resonance in a system with up to two critical points, and we describe a general technique that 
permits to treat other cases in a systematical manner. A separatrix divides the phase-space sections into regions of different topology: inside the 
libration region the evolutionary tracks bring the observed twin-peak frequencies to an exact rational ratio, whereas in the circulation region 
the observed frequencies remain off resonance. Our scheme predicts the power should cyclically be exchanged between the two oscillations. 
Likewise the high-frequency QPOs in neutron-star binaries, also in black-hole sources one expects, as a general property of the non-linear 
model, that slight detuning pushes the twin-peak frequencies out of sharp resonance. 
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1. Introduction 

Twin peaks occur in high-frequency (~ 50-450 Hz) power 
spectra of X-ray (~ 2-60 keV) lightcurves of several black-hole 
candidates (see van der Klis I2006I McClintock & Remillard 
2003 for recent reviews of observational properties and the- 
oretical interpretations). This transient phenomenon seems to 
be connected with the kilohertz quasi-periodic oscillations 
(QPOs) in neutron star sources, of which more examples are 
known (about the tens at present). The nature of black hole 
high-frequency QPOs remains puzzling despite variety of mod- 
els proposed in the literature. In neutron-star low-mass X-ray 
binaries these twin QPOs are known to occur often simultane- 
ously and they can be highly coherent (Q > 10 2 ; e.g. Barret 
et al. 2005 1, slowly wandering in frequencies between different 
observations, whereas in black-hole candidates the QPO co- 
herency appears to be lower (Q ~ 2-10) and the presence of a 
pair pops up only when a collection of observations is carefully 
analyzed. 

In Abramowicz & Kluzniak (20011 and Kluzniak & 
Abramowicz (2001 1 an idea of accretion disc resonance was 
proposed, which naturally incorporates pairs of frequencies oc- 
curring in a ratio of small integer numbers. This scheme pre- 
dicts the observed frequency ratios in black-hole QPO sources 



should prefer the 3:2 ratio; it also suggests this could be 
understood if a non-linear coupling mechanism operates in 
a black-hole accretion disc, where strong-gravity effects are 
essential. Indeed, especially in those black-hole candidates 
where the high-frequency QPOs have been reported, they oc- 
cur very close to the ratio of small integer numbers, 3:2 in 
particular (Miller et al. 120011 Strohmaver l2001l Homan et al. 
2004 Remillard et al. 120051 Maccarone & Schnittman|2005 1. 
Nowadays, the original account can be viewed as a naive 
model with the internal resonance. Various realizations of this 
scheme have been examined in terms of accretion disc/torus 



oscillations (e.g. Abramowicz et al. 2003 Bursa et al. 2004 



Kato 2004, Li & Naravan 120041 Schnittman & Rezzolla 2005 
Zanotti et al. l2005> . 

So far the "right" model has not yet been identified. 
However, it has been recognized that fruitful knowledge about 
common properties of high-frequency QPOs can be gained by 
investigating a very general resonance scheme, which likely 
governs matter near a compact accreting body. To this aim, 
Abramowicz et al. (2003) examined the epicyclic resonances 



in a nearly-geodesic motion in strong gravity. Rebusco (2004i 
and Horak ( 2004 1, by employing the method of multiple scales 
(Nayfeh & Mook 1979 1, have demonstrated that the 3:2 res- 
onance is indeed the most prominent one near horizon of a 
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central black hole. Only certain frequency combinations are 
allowed, depending on symmetries which the system exhibits, 
and only some of the allowed combinations have chance to give 
rise to a strong resonance. 

In the present paper we pursue this approach further and we 
find tracks that an axially symmetric system with two degrees 
of freedom, near resonance, should follow in the plane of en- 
ergy (of the oscillations) versus radius (where the oscillations 
take place). We show different topologies of the phase space in 
the way that closely resembles the method of disturbing func- 
tion, familiar from the studies of the evolution of mean orbital 
elements in celestial mechanics (Kozai 1962] Lidov|T962 1. The 
analogy is very illuminating and it provides a systematic way 
of distinguishing topologically different states of the system. In 
particular, one can discriminate regions of phase space where 
the observed frequency ratio fluctuates around an exact ratio- 
nal number from those regions where this ratio remains always 
outside the resonance. Our model suggests that even black-hole 
twin QPOs should vary in frequency and they should not stay 
at a firmly fixed frequency ratio, albeit the expected variation 
is very small - certainly less than what has been frequently re- 
ported in neutron star binaries and what can be tested with the 
data available at present. 

2. A conservative system with two degrees of 
freedom 

2.1. Non-linear terms in the governing equations 

Let us consider an oscillatory system with two degrees of free- 
dom, which is described by coupled differential equations of 
the form 



6p + cof. 6p - f p (6p, 60, dp, 60), 
60 + a)] 60 = f e (6p, 60, 6p, 59). 



(1) 
(2) 



We assume that the right-hand-side functions are nonlinear 
(their Taylor expansions start with the second order), and that 
they are invariant under reflection of time. Clearly, these equa- 
tions include the case of a nearly circular motion under the in- 
fluence of a perturbing force: 6p and 69 are small perturbations 
of the position, whereas co r (r) and u>g(r) have a meaning of the 
radial and the vertical epicyclic frequencies along the circular 
orbit r = ro, = 7r/2: 

d 2f U , 1 ( ,. 



(4) 



r dr 2 ' r 2 \ dO 2 
where the effective potential is 

I 2 

tl(r, 0) = <D(r, 9) + 



2r 2 sin 2 6>' 

<I>(r, 9) is the gravitational potential, for which axial symme- 
try and staticity will be assumed. These assumptions make 
our system qualitatively different from models requiring non- 
axisymmetric perturbations. 

As a generic example, let the total gravitational field be 
given as a superposition of the central potential of a spherical 
star, €> s (r), plus an axisymmetric term <D r (r, 9), 



For the central field we assume the form 
GM 

O s (r) = - — , 
r 



(6) 



where we set f — r or T — r - R s (to adopt the Newtonian 
or the pseudo-Newtonian approximations; R s = 2GM/c 2 ). We 
assume a circular ring (mass m, radius a) as a source of the 
perturbing potential, 



*r(r,fl) = - 



2Gm K(k) 



(7) 



where K(k) is a complete elliptical integral of the first kind, 
B(r, 0) = r 2 +a 2 + 2ar sin 0, k{r, 0) =4ar sin 0/B(r). 

At this point a remark is necessary regarding the interpre- 
tation of the potential l|5}-0: we conceive it as a toy-model 
for strong-gravity effects and internal resonances that we seek 
in the system of a black hole and an accretion disc, but the ori- 
gin of the perturbing potential €> r (r, 0) is not supposed to be 
the gravitational field of the accretion disc itself. Of course the 
inner disc is not self-gravitating in black hole binaries and the 
ring is not introduced here with the aim of representing the ac- 
cretion flow gravity. What we imagine is that hydrodynamic 
and magnetic forces are producing qualitative effects, which 
can be captured by the ring potential in our equations Q-Q. 
On the other hand, our approach is rather general and it is worth 
remembering that the same formalism can be successfully ap- 
plied also in systems where the disc self-gravity plays a non- 
negligible role. 

Also in general relativity, weakly perturbed (i.e. nearly- 
geodesic) motion of gas elements orbiting around a 
Schwarzschild black hole can be described by effective poten- 
tial that contains a spherical term arising from the gravitational 
field of the central black hole, plus a perturbing term, which we 
assume axially symmetric. Naturally, the interpretation of the 
perturbing term is more complicated if one would like to derive 
its particular form from an exact solution of Einstein's equa- 
tions. We do not want to enter into complications of a special 
model here (but see Letelier 2003 , Karas et al. 2004, Semerak 
2004 and references cited therein for a review and examples of 
spacetimes that contain a black hole and a gravitating ring in 
general relativity). Also, one may ask whether the adopted ap- 
proach can comprehend frame-dragging effects of Kerr space- 
time while the black hole rotation is considered as perturbation 
to the spherical field of a static non-rotating black hole; even 
if a general answer to that question would be affirmative, the 
Kerr metric has rather special properties concerning the inte- 
grability of the geodesic motion and the form of non-sphericity, 
which quickly decays with radius. Our current opinion is such 
that the rotation-related effects of Kerr metric cannot be viewed 
as the origin of the perturbation required for black-hole high- 
frequency twin QPOs. 

The angular momentum of a test particle orbiting the centre 
on an equatorial circular orbit is 



fer 3 



d<S> , 

-5- I = GMr 

or 



r z (r + a)E(ko) + (r - a)K(k ) 



n(r 2 — a 2 ) 



(8) 



<D(r,0) = <D s (r) + $ r (r, 0). 



(5) 



where E(k) is a complete elliptical integrals of the second kind 
and the right-hand-side terms are evaluated at the orbit radius, 
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Fig. 1. Angular momentum £{r) (thick line) of a test particle on a circular orbit in the combined gravitational field of a spherical 
body and a ring. The mass of the ring relative to the central mass is /j = 0.1 and the particle mass is set to unity. Left panel: the case 
of Newtonian central field. Right panel: the pseudo-Newtonian case. Radius has been scaled with respect to the ring radius (here, 
a = 9/?s), and the angular momentum has been scaled by the value of Keplerian angular momentum ^ic(a)- Keplerian angular 
momentum in the central field is also plotted (thin line). Circular orbits are Rayleigh unstable and the epicyclic approximation is 
inadequate in regions where the angular momentum decreases with radius. 




Fig. 2. The behaviour of epicyclic frequencies (thick lines; u>e - top, u r - bottom) provides basic clues to the origin of oscillations 
and the effects expected to occur in strong gravity. Left/right panels refer to the Newtonian/pseudo-Newtonian cases, as in the 
previous figure. Units of the ordinates are scaled by ^GM/a 3 to make frequencies dimensionless. The frequencies in absence of 
the ring (/i = 0) are also indicated for reference (thin lines; u>g(r) — a> r (r) in the Newtonian case). 



ko=k(ro, tt/2), fism/M. Figure^captures a typical curve of the 
angular momentum as a function of radius. The corresponding 
orbital frequency is 



2 (r) = — 



(r + a)E{k {) ) + (r - a)K(k ) 



n(r 2 — a 2 ) 



Equation gives the epicyclic frequencies, 



2, v GM 
w r W = — r 

2 GM 
0Jg(r) = — 



r 2 (3f - 2r) 2/j. (r - a) 2 K(k ) - a 2 E(k ) 
r 5 n r 2 (r — a) 2 (r + a) 



r 1 + 2pL 



E(k ) 



7i (r - a) 2 (r + a) 



(9) 

,(10) 
(ID 



The epicyclic frequencies as functions of radius are plotted in 
Figure |2 (we will drop the index "0" for the sake of brevity). 



The difference between the radial epicyclic frequency and the 
orbital frequency gives the shift of pericentre. The difference of 
vertical epicyclic and orbital frequencies gives the nodal pre- 
cession. 

Internal resonances can occur in the system Q-(|2|i. In or- 
der to capture this phenomenon, we carry out a multiple-scales 
expansion (Nayfeh & Mook 1979) in the form 



Sp(t, e) = Y J fpnVJ, m, e) = 2 fOnVJ, 



(12) 



n=\ 



where T M = e^t are treated as independent time scales. We will 
terminate the expansion at the fourth order (ji = 0, 1, 2, 3; the 
number of time scales is the same as the order at which the 
expansions are truncated). 
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Time derivatives take a form of expansions 



4 4 



H=0 v=0 



p=0 



(13) 



where Dp- = d/dTp. The method tackles the governing equa- 
tions in their general form. 

As a specific example we can adopt an explicit form de- 
scribing the orbital motion, 



6p + (o 2 r 6p = (l + 6p)69 2 - 
„ 5p66 



60 + 0)1 



1 +6p 



r dr 
1 



- ai r dp 



(l+Sp) 2 86 



ul 60 



(14) 
(15) 



Because equations (I14i — d!5t are conservative, the growth of 
energy in one mode must be balanced by the energy loss in the 
other mode. Close to resonance radii, where the two epicyclic 
frequencies are in ratio of small integers, the periodic exchange 
of energy should occur in a more pronounced rate. Because am- 
plitudes of the oscillations are connected with eccentricity and 
inclination, the solution alternates between an inclined, almost 
circular trajectory at certain stages, and an eccentric, almost 
equatorial case at some other time. We remark that, in a non- 
linear system, the eigenfrequencies u> r , a>g are expected to dif- 
fer from the observed (i.e. predicted) frequencies, which can be 
revealed e.g. by Fourier analysis of data time series. Relevance 
of this fact for QPOs was first recognized by Abramowicz et 
al. (|3003l and Rebusco (2004), when they discussed a model 
for Sco X-l. It was further employed by Horak et al. (|2004i, 
who suggested that a connection should exist between the high- 
frequency QPOs and normal-branch oscillations. 

We expand the effective potential derivatives and the func- 
tions f p and fe into Taylor series, up to the fourth order about 
a circular orbit. The expansion provides many nonlinear terms 
containing various derivatives 



MJ ~WdeJ) M2] 



(16) 



By imposing constraints on the potential and its derivatives, 
we identify resonances that are expected in a particular system 
(and we reject those resonances that cannot be realized). In the 
next subsection, formulation of these constraints is still kept 
completely general, valid for the arbitrary form of H. Only 
later, in subsection 12.31 we come to our original motivation 
from the orbital motion around a black hole and we employ the 
symmetry of the potential U. We consider the case of potential 
symmetric with respect to the equatorial plane. This implies the 
condition w,(2i+i) = 0, k e N, somewhat reducing the number 
of terms in the expansions. 

Amplitudes of the oscillations are characterized by a small 
parameter: bp ~ e, 69 ~ e. We impose the solvability con- 
straints and seek a solution in the form dl2> . To this aim, we 
first write the explicit form of these constraints in different or- 
ders of approximation. 



Table 1. Possible resonances and secular terms in the second 
order of approximation. Only regular secular terms are present 
if cog and w, are outside resonance, as indicated in the left col- 
umn. Subsequent lines refer to a system in 1:2 and 2:1 reso- 
nances, respectively. In each record corresponding to the res- 
onances, the first/second row gives an expression for secular 
terms in radial/vertical oscillations. To simplify our notation 



we introduced A a = cf 1 and K a = C^' . 



co e :co r 


Secular terms 


Outside 


-2icn r DiA p 




resonance 


-2lWgD[Ag 




1:2 


—2\u),-D\A p , 


KiootA-pAg 




— 1\(i)qDiAq, 




2:1 


—2\u},-D\A p , 


Koo02^ 




-TXlDgDyAg, 


A-onoApA-g 



2.2. Solvability constraints 

In the first order, we obtain equations describing two indepen- 
dent harmonic oscillators, 



(d 2 + w 2 ) Pi =0, {Dl + (4)9 X =0. 
The solutions can be expressed in the form 
pi = A p + A_ p , 9i = A g + A_ fl , 



(17) 



(18) 



where we denote A x = A x e l "" T °, A- x = A'e -1 "^ 7 * with x = p 
or 9. The complex amplitudes A x generally depend on slower 
scales, T\, T2 and T3. 

An algorithmic nature of the multiple-scales method allows 
us to determine the form of solvability conditions in a conserva- 
tive system with two degrees of freedom. The conditions arise 
by eliminating the terms that are secular in the fastest variable, 
To. This constraint is imposed, because otherwise the solutions 
in the form of series would not converge uniformly. In fact, the 
reason why the same number of scales is required as the or- 
der of approximation in the expansion (I12> is that one secular 
term gets eliminated in each order, and therefore the functions 
A x {Tn) are determined by the same number of solvability con- 
ditions as the number of their variables. 

In the second order, the terms proportional to e 2 in the ex- 
panded left-hand side of the governing equations are 

[S'x + cj 2 x x] 2 = (pi + co 2 .) x 2 + 2ico x DiA x - 2ioj x D x A- x . (19) 

On the right-hand side, the second-order terms result from the 
expansion of the nonlinearity f x {5p, 66, Sp, 69) with A)Pi an d 
Dq6\ in place of 6p and 66, respectively. They can be expressed 
as linear combinations of quadratic terms constructed from A ±p 
and A ± g\ 

[f x (6p, 69, 6p, 69)] 2 = J] C^A%A^A a _^ , (20) 

\a\=2 
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where a = (a\, . . . , 0-4) and \a\ = a\ + . . . + a^. The constants 
are given by angular frequencies of at x and by coefficients 
of the Taylor expansion of f x . Equating the right-hand sides of 
eqs. (THfl-CT we find 



+ City X 2 



(21) 



-2ico x D\A x + 2\(D X D\A- X 

, V r^T" 1 ~A a2 A hi A h4 

The rhs of eq. (12 1 1 contains one secular term independently 
of the eigenfrequencies a> r and a>g. However, additional sec- 
ular terms may appear in the resonance. For example, when 
10,- * 2u>e, the terms proportional to A 2 , in the p-equation 
(x — » p) and A p A_ g in the 0-equation (x — > 6) become secular 
and they should be also included to the solvability conditions. 
The similar situation happens when <±> r as a>g/2. These are inter- 
nal resonances, which show a qualitatively different behaviour: 
the corresponding terms are secular only for special (resonant) 
combinations of io r and ajg, contrary to the terms that appear 
always and are referred to as regular secular terms. Possible 
resonances in the second order of approximation and the cor- 
responding secular terms in eq. ( 12 11 1 are listed in Tabled Let 
us assume, for a moment, that the system is far from any reso- 
nance. Then 



D\A X = 0. 



(22) 



The frequencies and amplitudes are constant and the behaviour 
of the system is almost identical as what one finds in the linear 
approximation. The only difference is the presence of higher 
harmonics oscillating with frequencies 2ai r , 2a>g and \u) r ± u)g\. 
They are given by a particular solution of equation ( 12 U after 
eliminating the secular term, 



*2 = X q£a%a?a%a?, 

M=2 



(23) 



Under the assumption of time-reflection invariance, constants 
Qa are real and their relation to constants C„ becomes obvi- 
ous by substituting X2 into equation ( 12 li . We find 



(x) 



c 



w 

klmn 



[(k — l)a) r + (tn — n)u)e\ 



(24) 



In the third order, the discussion is analogous in many re- 
spects. The terms proportional to e 3 , which appear on the lhs 
of the governing equations, are given by 

\6x + co 2 x xj = {pi + cufj x 3 + 2\u x D 2 A x - 2\oj x D z A- x . (25) 

The terms containing D\X\ and D\x 2 vanish in consequence 
of the solvability condition \22\ . The rhs contains now cubic 
terms of the Taylor expansion. We obtain 

(Dq + w^.)x 3 = -2\cd x D 2 A x + 2\(d x D 2 A- x 
Zj a -p p A -e A e > 



(26) 



\a\=3 



where constants C„ are real. The secular terms are summa- 
rized in Table |2 together with resonances possible in the third 
order of approximation. Again, far from any resonance we 



Table 2. Possible resonances and secular terms in the third or- 
der of approximation. Individual records have a similar mean- 
ing as in Tab. 1 . 



dJg'.OJr 


Secular terms 






Outside 


2uo r D 2 A p , K 12 ookM a p 




, Kom|A 2 |A p 




resonance 


2\u)„D 2 A g , A 110 i|A p | 2 A fl . 


, A oi2|A 2 |A(, 




1:3 


2ioj r D 2 A p , K 120 o|A p | 2 A p 


, Koui|A 2 |A p , 


Koo3oA^ 




2iwf>D 2 Af>, A 110 i|A p | 2 A fl . 


, A oi2|A||A(j, 


A01 20 A p A_ g 


1:1 


2iw,.Z) 2 A p ^K 1 2oo|A p | 2 A e 
Kooi2|A e | 2 A fl , K 2ioA^,A_ 


, KomlAjlAp, 
e, K[ o2A_ p A e 


Ki U o|A p | 2 A fl , 




2ioj D 2 A e , A I10 i|A p | 2 A s . 


, A oi2|A 2 |A(j, 


A 2 ioo|A p | 2 A p , 




A o2i|A fl | 2 A fl , A1002A_ p 


A 2 , AoaioApA; 


-e 


3:1 


2icn r D 2 A p , Ki2oo|A p | 2 A p 


, Kom|A 2 |A p , 


K2ooiA 2 p A s 




2itD 9 D 2 A 9 , A U oi|A p | 2 A fl . 


, A 00 i2|A 2 |A e , 


Ao30oA^, 



eliminate the terms that are secular independently of a> r , a>g. 
The resulting solvability conditions take the form 

D 2 A P = [K 1200 |A p | 2 + K 0111 |A 9 | 2 ]a p , (27) 

D 2 Ag = [Aik,i|A p | 2 + A 00 i2|A fl | 2 ]A fl . (28) 

A particular solution of eq. \26\ is given by a linear combina- 
tion of cubic terms constructed from A ±p and A±g, 



*3 = £ Qa X) A%A p -A%A a e \ 

where all coefficients Q ( „' A) are now real. 

Finally, in the fourth order of the approximation, 

\8x + w 2 xj 4 = (pi + w 2 ) X3 + 2D 3 D xi + 2D D 2 x 2 . 



(29) 



(30) 



The operator DoD 2 acts on x 2 , given by eq. d23l >. The resulting 
form is found by employing the solvability conditions J27b and 
S3: 



2D D 2 x 2 = oj] £ J^A^-A^ 



(31) 



M=4 



where j£' are real constants. By expanding the rhs we arrive at 
the governing equation 



[d 2 + oj x ] 



-2\u) x Dt,A x + 2iw v Z} 3 A_ x 
Zj a -p p -e^s > 



(32) 



|tr|=4 



with C a real constants. Only one secular term independent of 
a> r and cog appears on the rhs: -2\(jj x D^A x . The sum contains 
only terms that become secular near a resonance. These terms 
and the solvability conditions are listed in Table|3] 
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Table 3. Possible resonances in the fourth order of approxima- 
tion. 



w g :u) r 


Secular terms 


Outside 


2ito r D 3 A p 




resonance 


2ico e D 3 A e 




1:4 


2ia> r D 3 A p , 
2ia> e D 3 A e , 




2:3 


2itL>,D 3 A p , 


K-mioApA^g 




2lOJgD i Ag, 




3:2 


2iu) T D 3 A p , 

2\lligDiAg, 


K2002^ p ^ 


4:1 


2iaj r D 3 A p , 


— 3 "a~ 

Ko3oiA p A s 




2l(jJgDiAg, 


Ao400^p 



A notable feature of internal resonances k:l is that ka> r and 
lu)e need not be infinitesimally close to each other, as might 
be expected from the linear analysis. Consider, for example, 
an internal resonance 1:2, i.e. a>g m 2a> r . By eliminating the 
secular terms, we obtain solvability conditions (see Tab.^ 



2\u) r DiA p + KiooiA_ p A e = 0, 

-HlOgDiAg + A()200A^ = 0. 



(33) 
(34) 



In each of these equations the first term is regular, while the 
second term is nearly secular (resonant) one. The solvability 
conditions give us the long-term behaviour of the amplitudes 
and phases of oscillations. Suppose now that the system departs 
from the sharp ratio by small (first-order) deviations u>g = 2a> r + 
ecr, where cr is the detuning parameter. The terms proportional 
to A-pAg and A p still remain secular with respect to the variable 
Tq. This can be demonstrated from A- p Ag = A* p Age 
A*Age [crT ' e [oJl To . Similar relation holds forA^. 



i(a> e -b> r )T 



2.3. Solution of the solvability constraints 

By comparing the coefficients with same powers of e on both 
sides of the Taylor-expanded governing equations we 
obtain relations for functions Pi(Tj) and 6i(Tj) that can be 
solved successively. After rearranging to a "canonical" form, 



[D 2 + oj 2 ]p n = ^ Kjj k iA'_ p A J p A k _ g A l e , 
[Dl + u> 2 g ]e n = ^ Ai jk ifi_ f fi p A k _ g A l g, 



(35) 
(36) 



where n is the order of approximation. In this way we iden- 
tify constants Kim and Ay^. The studied gravitational poten- 
tial is symmetric with respect to the equatorial plane, therefore, 
the series d35l— d36l cannot contain terms proportional to odd 
derivatives of the effective potential with respect to 9. Hence, 
contrary to a general case, only specific resonances occur here: 



a>g\u r — 1:2, 1:1, 3:2 and 1:4. These are all possible combi- 
nations that may occur within the given order of approxima- 
tion (the first three cases were originally identified by Rebusco 
2004 although she does not mention the fourth possible com- 
bination). A general argument of non-linear analysis suggests 
that the dominant resonances are those ones which correspond 
to ratios of small natural numbers, although not every con- 
ceivable resonant combination comes up in a given physical 
system. Indeed, it appears that the 3:2 ratio is the most im- 
portant case when the high-frequency QPO pairs are debated, 
however, the true role of this resonance has not yet been un- 
derstood; see also the discussion in Bursa (2005) and Lasota 
(2005). Proceeding further to higher-order terms of the expan- 
sions reveals even more resonances, but these are expected to 
be very weak (recently various kinds of weird combinations 



have been examined by Torok et al. 2005). Hereafter we con- 
centrate on the first three combinations. 



The case of 1 :2 resonance 

The solvability conditions take the form (cp. Tab. [l] and Horak 

I2004T) 



D\A p - —- — Kooo2^fl, 
D\Ag = -— !— AoiioApA_ e , . 

2tOg 

The coefficients of the resonant terms are given by 

rr _ 2 "12 A 2 M 12 

JM)002 - -Wg ~ ~ — > A 0110 - -^g > 

2r r 



(37) 
(38) 

(39) 



satisfying a mutual relation 2Kooo2 = Aono- 

The case of 1 :1 resonance 

The solvability conditions for the first order, D\A p - D\ Ag = 0, 
imply that the complex amplitudes A p and Ag depend only on 
the second scale T 2 . The 1:1 (a> r « a)g) resonance is the only 
epicyclic resonance of the system with reflection symmetry, 
which occurs in the third order of approximation. The depen- 
dence on T 2 implies a slower behaviour. The solvability condi- 
tions are 



D 2 A P = - 



2co r 



K1200 \Ap\ A p 



+K0111 \A g \ Ap + K 10 o2^-p A fl 



D 2 Ag = - — 
2a>g 



A 



1101 |A p | Ag 



+A0012 \Ag\ 2 A g + A 02W A 2 p A_ 6 



and the coefficients of the resonant terms are given by 



K 



1200 



0111 



2 M 40 

6u>g 



2u\ 2 



(40) 



(41) 



(42) 
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-6r u 30 + u l2 



8 3 M30 

- + — r 



K1002 = 7 
o 



-6i0n + 



6m 12 n r, «12«30 X 

- 3 M22 - 2 r u 3 o 5— 



2 2 



w; 



5 k? 



10 



«04 7mj2 -"*12 2 

Ao012 = -55 " + 6^1 + t w « 



A02IO - Ki002, 

Anoi = Koin. 

The case of 3:2 resonance 



(43) 
(44) 

(45) 

(46) 
(47) 



The solvability conditions involve both the third and the fourth 
orders. Hence, the amplitudes A p , Ag are functions of both time 
scales Tt, and T4. The elimination of regular secular terms in 
the third order (3u> r w 2u>g; see Tab.|2ji gives 



1 

2a> r 
i 

2u( 



DiA p = - — 
D 2 A g = - 
with the coefficients 

K1200 



K1200 A p A p + Koin \A g \ 2 A p 



Anoi |Ap| 2 A e + Aooi2|A fl | 2 A 6 



; I ^ M 30 1 ^ 

f"n I — ^«40 



K0111 - — 



1 we 

2 9 M12 
-15 w fl + — — — - - 4 M22 



(48) 
(49) 

(50) 



■18r M30 + 9Mi2 



n " ,2 



A0012 
Anoi 



M 4 135 u\ 2 _ 153 Ml2 135 2 
'27! + "64 ^2 ~ "lT + ~16 



'0 
Koin- 



(51) 

(52) 
(53) 



The elimination of the resonant terms gives the solvability con- 
dition in the fourth order (Tab. [5}, 



£>3 Ap = - K2002 A 2 ^ , EhA e 
where the resonant coefficients are 



A 3ioA!A-e, (54) 

2(jDg P 



,.2 ,.3 

15 . .2 , 27 ffU ^ 135 "l2 _ 243 M 12 

16 W fl + 32 



'0 



64 2 2 128 3 4 



"R "22 + Tfi 



27 M 12 "22 27 



16 



— nmo + 



M 12 M30 



64 



_9_ r M 2 2 «30 _ _8l_ r Q M 30 2 _ _gj_ r «12«3o 
16 



W7, 



256 ,2 



512 4 



1 92 9 r 0Ml2 "40 
"4 r «30 - 64 ^"40 - T28 2 ' 



Ao310 - -|K2002- 

By introducing the detuning parameter, 
<te3 2 = ^cr 2 + e o" 3 , 



(55) 
(56) 

(57) 



0.12 



0.08 



0.04 



8 = const 




0.04 



0.08 



0.12 



0.16 



Fig. 3. Comparison between the analytical constraint 
S(a p ,ag) = const (an ellipse), derived in multiple-scales 
approximation, and the corresponding exact (numerical) 
solution. The curly curve is the numerical solution of the 
oscillation amplitudes a p (t), ag(f). The agreement between the 
two curves demonstrates that accuracy of the approximation is 
satisfactory over the entire time span. See the text for details. 



the solvability conditions adopt the explicit form 

2iw,D 3 A p = K 2002 (A* p f A 2 g e^ + ^\ 
2\a>eDiAg = 



A A^ A* „i(o H 2 7~2+o~i7"i) 

Ao3l()A p A fl e 



2\u r D 2 A p = 
2ia)gD 2 A g = 



K1200 |A p | + K iii \A g [ 



Ai 101 Ap +Aooi2|A e r 



Ap, 

Ag. 



(58) 
(59) 
(60) 

(61) 



Finally, by substituting a polar form of complex amplitudes, 
A p = iape 1 ^ and Ag = \age l ^ b , we get a set of eight equations 
governing the long-term behaviour of phases and amplitudes: 



D 2 a p = 0, D 2 ag = 0, 

„ . K2002-2-2 • 
D 3 a P = TT7~ a o a 6 sm 7 



D^ag 



I6u> r p 6 " 16dJg 
D 2 (p p = [Knooflp + Koinflel > 

OU) r J 

D 2 <pg = -77— [Anoiflp + Aooi2fle] , 



A0310.3- 

a p a g sin 7, 



8 we 

K2002 2 
D3<P P = -— — a p a g cosy 



Di4>g 



A0310 -3 

■ a p cos y, 



(62) 
(63) 

(64) 
(65) 



16w r 16we 
where the phase function has been introduced as y(T 2 , r 3 ) = 
-cr 2 T 2 - 0-373 - 3(p p + 2<pg. The amplitudes a p and ag of the 
oscillations vary slowly, because they depend only on the third 
time-scale T3. Phases <p r and <pg of the oscillations evolve on 
both time scales T 2 and T3. 

3. The system evolution near the 3:2 resonance 

3.1. The integrals of motion 

The case of 3:2 resonance is particularly relevant for the high- 
frequency QPOs, both on observational and theoretical grounds 
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Fig. 4. The functions F\ = F(u), F2 = —F(u), and G(u) from 
eq. ( 18 It . The system evolution is limited to the interval {u.\,uq), 
where the condition \F(u)\ > \G(u)\ is satisfied. 



(cf. Abramowicz & Kluzniak 2001 Kluzniak et al. l2004l for 
arguments in favour of 3:2 ratio in high-frequency QPOs and 
for further references). We therefore discuss this case in more 
detail, however, similar discussion could be presented also for 
other resonances (Horak 2005i. Reintroducing single physi- 
cal time f, the equations for the second and the third order 
can be combined together. Time derivatives are then given by 
d/dt = e 2 £>2 + e 3 £>3. Amplitudes and phases of the oscillations 
are governed by equations 



a r - -^p<jj r a 2 .a 2 g s\ny, 

1 o, . „3 



--^Pojga'ag sin 7, 

tOg 

-o-ajg + — 



H r a 2 r 



2 a r I 

■ ^edg + — [c 



cosy 



(66) 
(67) 

(68) 



where a p - 



ea p , ag - eag and fi r , fig and /3 are defined by 
relations A 03 io = -5K2002 = fa 2 , K1200 - Anoi = o£jU r , and 
Kmi - A0012 = The amplitudes and phases are not mutu- 
ally independent; the equations <66>-(l67t imply that the quan- 
tity 



£ = a 2 + jal = const 



(69) 



remains conserved during the system evolution. Clearly, £ is 
proportional to the total energy of the oscillations. The exis- 
tence of this integral is a general property of a conservative 
system. Naturally, it is not limited to the particular form of the 
perturbing potential 0, which we consider here, and it holds 
in Newtonian, pseudo-Newtonian as well as general-relativity 
versions of the equations of motion (the pseudo-Newtonian 
case was examined, in detail, by Abramowicz et al. 2003J and 
Horak 2004 1. One can watch the accuracy to which £ is con- 
served in order to verify the analytical approach against the ex- 
act numerical solution (we show such comparison in Figure 
for the pseudo-Newtonian case of Abramowicz et al. l2003> . 

Equations J66l >— J67l > can be merged in a single equation by 
introducing the following parameterization: 



0.02 



0.01 



- 


1 V 

j/^TypeB/ \ \ 


_ 




1/ TypeC \ \ 






' \ 

Type A \ 


1 



0.5 



0.55 



0.6 



0.65 



0.7 



r / a 



Fig. 5. The 3:2 inner resonances in the gravitational field of a 
spherical pseudo-Newtonian star and a ring (a = 9/?s, = 0.1). 
The regions of different phase-plane topology are identified in 
(r, £)-plane. Three types can be distinguished according to the 
number of critical points: A - no critical point (the system is far 
from resonance); B - one critical point; C - two critical points. 

Then the oscillations are described by two equations for £(f) 
and y(f), 

i=±Poj e f(l-f)a 3/2 smy, (71) 
y = -crujg + ±u> e &\iJ. r f 

+ > (l - f) + ^(3 - 5£ 2 ) £ 1/2 cos y] , (72) 
which satisfy the identity 

£dy-rd£ = 0. (73) 

Substituting for £ and y from eqs. f7H-f72t. eq. ( I73> implies 
that 

T = 8 (w 2 )<r + £[/^->(W 2 ) 2 

+jS£ 3/2 £ 3 (l -f) cos y (74) 

is a second integral of motion. For a given value of energy £, 
the system follows T = const curves. In the other words, pro- 
jection of the solution onto (y, £)-plane satisfies 



T{y, £) = const. 



(75) 



This allows us to construct two-dimensional phase-space sec- 
tions in which the system evolution takes place. 

3.2. Stationary points and the phase-plane topology 

Stationary points are given by the condition ^ = y = 0. 
According to eq. ilH . y-coordinate of these points satisfies 
siny = 0, and therefore y - kn with k being an integer. 
Substituting y — and cos y - + 1 in eq. 1721 . we find a cu- 
bic equation, 

(70) -4 ( r+[^ 2 + |^(l-^ 2 )]£+^(3-5^)£ 3/2 = 0. (76) 
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The solution determines ^-coordinate of the stationary points. 
In the case of small oscillations (£ <k 1), the solution can be 
approximated by keeping only terms up to the linear one in & 
in eq. (I76> . We obtain 



where T can be roughly approximated as 



2 _ 9d--fi B 

b 9 

-y"r - He 



(77) 



where cr = cr/S. The first correction to this solution is of the 
order of £^ 2 . Deviations between ^-coordinates of stationary 
points at odd and even multiples of n are of the same order. 

The solution J77I lies within the allowed range provided 
that jfi r § fr § ifig with the denominator D = - // e § 0, 
respectively. This can be expressed in terms of energy S: given 
the detuning parameter cr, stationary points appear in the (y, 
plane if the energy of oscillations satisfies 



9 — $ & $ 4 — for D % 0. 

He H>- 



(78) 



The examination of phase-plane topology near critical 
points leads to the equation 



(dj_ 



idj 



d% dy 
dy di 







(79) 



for eigenvalues A of the system of linearized equations (17 1 i and 
J72> . Evaluating the partial derivatives at the critical point and 
keeping only the terms of the lowest order in £ 3 ^ 2 , we obtain 



■5/2 



(80) 



Examining the sign of A 2 demonstrates that the critical points 
of central topology alternate with those of saddle topology. 

3.3. The time dependence 

The equation for g(f) can be derived by combining eqs. 17 H 
and 11741 , Eliminating cosy, we arrive at the relation 



Id? = F 2 (u) - G 2 (u), 



(81) 



where we introduced a new variable u(t) = £ 2 . The constant "K 
and functions F(u) and G(u) are defined by 



&< 2 Wo/3 



(82) 



F(u) = M 3/2 (l -u) (83) 

G(u) s Wn- " 8cr(1 " u) ~ ^ + ^ eS(1 " M)2 J ■ (84) 

The motion is allowed only for u 2 > 0, and so the condition 
+F(u) = G(u) gives us two turning points, u\ and u%, between 
which the evolutionary path oscillates. The functions +F(u) 
and G(u) are plotted in Figure|4] 

The period of energy exchange can be found by integrating 
eq. EJ: 



du 



y/F 2 (u)-G 2 (u)' 



(85) 



16tt 



-3/2 



(86) 



Notice that this time-scale is longer than the period of indi- 
vidual oscillations and it says how fast the system swaps itself 
between the radial and vertical oscillation modes. The simple 
estimate d86l is quite precise in most parts of the phase space, 
although it becomes inaccurate near stationary points, where 
the rate of energy exchange slows down. 

Finally, we illustrate our results with a simple case, which 
we already introduced at the beginning of the paper Csect. l2.lt : 
the gravitational field generated by a pseudo-Newtonian star 
and a narrow circular ring. The resonant condition <i>o/co r = 3/2 
is now fulfilled at three different radii, two of them lying be- 
tween the star and the ring, and the third one outside the ring. 
The resonances occurring at first two radii are called the in- 
ner resonances, whereas the latter is the outer resonance (not 
to be confused with 'internal resonance', which they are all). 
We restrict ourselves to the inner resonances, for which we 
find cr, fi r , fie and /3 as functions of r. For a fixed radius, in- 
equalities J78i give us the energy range of the oscillations. The 
result is shown in Figure |5J where we identify three different 
phase-plane topologies in the (r, £)-section. These can be dis- 
tinguished by the number of critical points and the shape of 
separatrices. The topology change is evident in Figure[6] where 
two-dimensional plots are constructed for the integral of mo- 
tion T(£,y). 

3.4. Frequencies of the resonant oscillations 

Equations d63l >— d65l > give the shift of actual (observed) frequen- 
cies of oscillations, w* and of g , with respect to the eigenfre- 
quencies <±> r and ojg: 



OJ r — <D r + (p r , 0> g — (i)g + <pg. 

These relations can be combined to find 
2co* g - 3oj* — 2ioe - 3o> r + 24>e - 30 r = 7- 



(87) 



(88) 



The observed frequencies are in exact 3:2 ratio if (and only if) 
the time-derivative of the phase function y vanishes. An imme- 
diate implication for the frequencies of stationary oscillations 
with constant amplitudes is that they stay in exact 3:2 ratio, 
even if the eigenfrequencies may depart from it. Outside sta- 
tionary points, it is evident from Fig. |6] that y — represents 
turning points on libration tracks (those ones, which are encir- 
cled by the separatrix curve). Hence, eq. d88l discriminates be- 
tween librating and circulating trajectories in the (y, £)-plane. 
Circulating trajectories span the full range of —n < y < n and 
they do not contain any turning point; y remains nonzero and 
the twin frequencies never cross the exact 3:2 ratio in the region 
of circulation. On the other hand, there are two points y — on 
each librating trajectory. In such state the ratio of observed fre- 
quencies slowly fluctuates about 3:2. 
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Fig. 6. Different topologies of the phase-space sections are shown in (£, y)-plane. The system evolution follows the contour lines 
T(^, j) = const. Their shape determines a possible range of the frequency change of the QPOs. Separatrix curve (thick line; cases 
B and C) divides the regions of circulating tracks from the regions of libration. Parameters and notation as in previous figure: 
A (left; r = 0.6a, & = 0.02), B (middle; r = 0.65a, & = 0.02), C (right; r = 0.55a, 8 = 0.05). 



4. Conclusions 

We have discussed the resonance scheme for high-frequency 
QPOs via multiple-scales analysis, assuming an axisymmet- 
ric conservative system with two degrees of freedom. This ap- 
proach provides a useful insight into general properties that are 
common to different conceivable mechanisms driving the os- 
cillations, although it does not address the question how the 
observed signal is actually formed and modulated. In our sce- 
nario, amplitudes and phases of the oscillations are mutually 
connected and they follow tracks in the phase space with dis- 
tinct topologies. The particular form was assumed to couple 
the oscillation modes via the non-spherical terms in the gravi- 
tational field of a ring. We consider this to be a toy-model for 
rather general behaviour, which should take place in any sys- 
tem governed by equations of type Q-(|2j. 

We assumed the Newtonian (or the pseudo-Newtonian) de- 
scription of the central gravitational field with a perturbation 
by an aligned ring as an example. The adopted form is not es- 
sential for general conclusions. In fact, equations cover 
also the nearly-geodesic motion around a Schwarzschild black 
hole. Compared with the pseudo-Newtonian case, general rel- 
ativity does not bring qualitatively new features, as long as 
the system is conservative (additional terms will arise in the 
expansions, which then translate to slightly different value of 
the resonance radius and to different duration of time inter- 
vals in physical units). A natural question arises whether the 
gravitational field of a rotating black hole could provide the 
perturbation required for the internal resonance in a surround- 
ing disc. We considered this possibility, however, it is unlikely 
that Kerr metric could by itself suffice: in the weak-field limit 
non-spherical terms seem to be incapable of creating separa- 
trices in the phase-space sections discussed above, whereas in 
the full (exact, vacuum) Kerr metric the special mathemati- 
cal properties of the spacetime ensure the integrability of the 
geodesic motion, and hence prevent the occurrence of inter- 
nal resonances. Therefore, the problem of a specific mecha- 
nism launching and maintaining the oscillations remains unan- 
swered. 

Various options for the generalisation of our scheme could 
be motivated by papers of other authors who proposed specific 



models including non-gravitational forces (see Petri 2005 for a 
recent exposition of the problem and for references). As a next 
step towards an astrophysically realistic scheme, one should 
take dissipative and non-potential forces into account, as well 
as non-axisymmetric perturbations. These will allow our sys- 
tem to migrate across contours in the phase-plane and to un- 
dergo transitions when crossing separatrices. Such additional 
terms could also supplement the influence of external forcing 
and initiate the oscillations of the system. The internal reso- 
nance would then define the actual frequencies that are excited; 
this way the strong gravity unmasks itself. 
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